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ABSTRACT 

We use an empirical approach to model the stellar mass of galaxies according to their host 
dark-matter haloes and subhaloes ('HASH' models), where each galaxy resides in a subhalo 
taken from a large iV-body cosmological simulation. This approach allows us to study the 
mass relation between subhaloes and galaxies (MR) using various observational constraints 
at redshift zero: the weak lensing signal (WL), the two point auto-correlation function (CF), 
and the stellar mass function of galaxies (SMF). Our method is based on modeling the lensing 
signal directly from the cosmological TV-body simulation, and should thus be more accurate 
(at least for massive objects) than other methods based on analytic halo models. We find that 
the WL does not provide a strong constraint on the MR. The current observational accuracy 
allows for more than a factor of 10 freedom in the subhalo mass of central galaxies, for a given 
stellar mass. The freedom for satellite galaxies is much larger, providing a very poor constraint 
on the number-fraction of satellite galaxies (0.05 - 0.8). These results are not significantly 
modified when using both the SMF and WL as constraints. We show that for the most massive 
galaxies, observational constraints based on the CF with 0. 1 dex errors, are equivalent to 0.05 
dex error in the WL. For intermediate and low mass galaxies the WL and CF constrain the 
MR in a different way. Although the WL is currently not adding much information at these 
masses, it has the potential of being important using future, more accurate measurements. 
The models found here do not match simultaneously the observed CF and WL signals, and 
show a limited ability to match the WL & SMF. We suspect that this is partially due to the 
cosmological model assumed here and we therefore adopt a mock WL signal through most 
of this work. In comparison to previous models in the literature, the method presented here is 
probably more general, as it reveals a larger range of solutions for a given set of observational 
constraints. 

Key words: galaxies: abundances; galaxies: formation; galaxies: haloes; galaxies: mass func- 
tion; galaxies: statistics; cosmology: large-scale structure of Universe; gravitational lensing: 
weak 



1 INTRODUCTION 

The relation between haloes and galaxies marks the crossroad of 
two major pathways. On the one hand dark-matter is accurately 
modeled by large TV-body cosmological simulati ons, where the 
properties of haloes are quantified in detail (e.g. ISpringel et al.l 
l2005t l iKlypin et alj|2011f) . Moreover, the large dynamical range in 
these simulations provides a full description of the dark-matter sub- 
structure that is believed t o host galaxies toiemand et al J [20071 : 
iBovlan-Kolchin et alj[200gh . On the other hand, the properties of 
luminous galaxies can now be measured over large cosmologi- 
cal volumes, redu cing the statistica l uncertainties to a le vel of a 
few per cent (e.g. iDrorv et alj|2009l ; iBernardi et alj|2010h . Obser- 
vational techniques for estimating physical properties of galaxies, 
like stellar mass and star-formation rate, are becoming more and 
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more reliable (e.g. lNordon et alj2012l) . making it feasible to model 
the physical quantities of galaxies directly. 

Various different methodologies have been developed in order 
to study the mass relation between galaxies and haloes (hereafter 
MR). Such models need to be flexible enough, so they can fit the 
observational constraints to a high accuracy. They also need to be 
computationally efficient, to allow quick exploration of the parame- 
ter space. Lastly, the model parameters should be directly related to 
the observational constraints (e.g. the stellar mass function and the 
auto-correlation function). Due to the above, all the models used 
to study the MR in detail are based on empirical relations between 
haloes and galaxies. However, in order to make reliable predictions 
on the MR, we need our models to be general enough, so they can 
encompass the full ensemble of solutions. One of the main targets 
of this work is to emphasize the importance of this issue. 

Halo occupation distribution models (HOD) were the first 
to make strong predictions on the MR, and are being widely 
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us ed ever since (e.g. Ijing et alj 1 19981 ; Coorav & Sheth |2002|; 



Berlind & Wein berg! 120021. iTinker et al. 120051: IZehavi et alJ|2005l : 
Zheng & Weinberg! |2007| ; ITinker et alj [20081 ; IZehavi et al.N201lh . 



These models assume that each halo hosts a number of satellite 
galaxies N g , where N g follows a power law as a function of the 
host halo mass. The galaxies are positioned at random locations 
within the halo, with a probability function that follows a scaled 
version of the dark-matter density profile. In order to predict ob- 
servables like the abundance and clustering of galaxies, HOD mod- 
els use the number density of haloes and their clustering properties 
by fitting results from large iV-body simulations. 

Since HOD models are based on analytical arguments, 
they are very useful for developing an insight on the model 
ingredients, and how the y affect the observed quantities (e.g 
iBerlind & Weinberg! 2002h . In addition, HOD models are relatively 
simple to understand, as they include just a few free parameters. 
These benefits also allow for detailed studies on ho w the underlying 
cosm o logical paramete rs affect the model r e sults {Abazaiian et al.1 
|2005| ; lYoo et alj 120061 ; ITinker et al.1 12007| ; iBaldauf et al.1 12010|) . 
which are then used for making predicti ons for the cosmological 
parameters (e.g. lMandelbaum et alj|2012h . 

A variant of HOD, based on the cond i tional luminosity func- 
tion (CLF), was developed by lYang et al] {2003). Within this ap- 
proach the properties of galaxies (e.g. the luminosity function) 
are matched for a given range of halo mass, using analytical fit- 
ting functions for the contribution from central and satellite galax- 
ies. Although this parametrization is somewhat different than the 
usual HOD models, and includes different degrees of freedom, the 
two approaches are very similar. CLF studies were used to in- 
terpret v arious types of observation s (e.g. mass-to-light ratios in 
clusters. Ivan den Bosch et aUl2003l), and to constrain the cosmo- 
logical parameters dvan den Bosch et al. | |2012| ; lMoreetalj|20T2l ; 
ICacciato et all2012h . 

Lastly, abundance matching models (ABM) populate galaxies 



mological simulations 1 Vale & Ostriker 2004; Conrov et alj 2006; 


Shankaretal. 2006; IConrov & Wechslei 


20091; Behroozi et al. 
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Guoetal. 2011; 
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20121; iReddick etal] 12012). In standard ABM. the MR between 



galaxies and their host subhaloes is set by matching the luminosity 
function (or stellar mass function) of galaxies to the mass function 
of subhaloes. The only important detail lies in the way the subhalo 
mass is defined. While the mass of central subhaloes is derived at 
the observed redshift, the mass of satellite subhaloes is set by their 
mass just before becoming a satellite, which typically occurs a few 
Gyr before the observed epoch. ABM models have no free param- 
eters, and seem to be too simplified at first sight. However, these 
models are usually very successful in matching the mass function of 
galaxies (by construction), their auto-corre lation function, and ev en 
weak gravitational lensing measurements iTasitsiomi et alj|2004l) . 

All the models above use some substantial assumptions that 
might restrict the nature of the accepted solutions. Both HOD and 
CLF are based on specific parametrization. For example: the aver- 
age number of galaxies within a given halo mass (N g ), the func- 
tional shape of the conditional luminosity function, and the lo- 
cation of satellite galaxies. Although most of these were tested 
against hydrodynamical and iV-body simulations or group catalogs 
(Rravtsov et al. 2004; Zheng et al ] |2005l ; iMandelbaum et alj|2005t 
ISimha et al.l2012h . these simulations do not necessarily span all the 
possible physical scenarios of galaxy formation. This problem is 
much more severe in ABM models, as they assume a specific rela- 



tion between satellite and centra l galaxies, forcing bot h populations 
to follow the same MR (see e.g. lNeistein et alj|2011al) . 

The concern that some assumptions are too restrictive also 
arises when com parin g different s t udies. For example, both 
lYang et alj J2012b and iMoster et al. I d2012l) have tried to self- 
consistently incorporate the evolution of the MR with redshift. Al- 
though these models assume a very different behavio ur for satellite 
galax ies (satellites are not allowed to form stars in IMoster et al.1 
I2012L but grow significantly in Yang et al. 2012) both models are 
able to fit the data. This example indicates that more complex mod- 
els do not necessarily capture all options, and some non-negligible 
freedom might be found after changing the assumptions of each 

model above. 

In our previous study dNeistein et alj|20~l lbl hereafter paper- 
/) we have proposed a synthesis between HOD and ABM termed 
'HASH' (an acronym for 'halo and subhalo'). Within this approach 
we assign galaxies to subhaloes that are taken from a large TV-body 
simulation. We allow the mass of satellite galaxies to depend on 
both the halo and subhalo masses. In this way we use a more reli- 
able estimate for the number and location of satellite galaxies than 
HOD, while allowing the halo mass to play a role in shaping the 
clustering of galaxies. In our implementation, we try to use the 
minimal set of assumptions possible: the mass relation for central 
galaxies is not restricted to follow any functional shape, the de- 
pendence on both the halo and subhalo masses is being constrained 
only slightly, the number of satellite subhaloes can change by using 
different dynamical friction estimates, and the location of satellite 
galaxies can be modified. All these ingredients aim at making our 
approach more general, and should allow us to find more possible 
models that fit the data. 

In paper-/ we have shown that our approach points to a 
large degeneracy in the models that fit both the stellar mass 
function of galaxies (SMF), and their projected two-point auto- 
correlation function (CF). The natural next step would be to 
check how additional constraints of very different nature can re- 
strict the degeneracy of our models. In this work we add mea- 
surements of weak gravitational lensing (WL) as it is claimed 
to provide str ong co nstrains on models when combined with the 



CF JYqo et alj|2006l: ICacciato et ai1l2009l: Eeauthaud et al 
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Unlike previous studies, we will show below that using observa- 
tional constraints from WL and SMF can hardly limit the mass re- 
lation between haloes and galaxies. Using all constraints together 
(i.e. WL, SMF and CF) is currently equivalent to using only the 
CF and SMF. However, we show that future WL observations with 
increased accuracy at small stellar masses should be able to signif- 
icantly contribute to our knowledge on the MR. 

The general principle of the HASH approach is motivated by 
the properties of galaxies within a sem i-analytic model (SAM). We 
have shown in iNeistein et alj j201 lal) that the combination of the 
subhalo and halo mass can reproduce the clustering of the SAM 
galaxies to a reasonable accuracy. In general, the subhalo mass is 
more closely related to the stellar mass, and gives a good handle 
on the location of galaxies, while the halo mass strongly affects the 
clustering prop erties of galaxies, and induc e possible e nvironmen- 
tal eff ects (e.g. [K hochfar & Ostrikeill2008l) . Recently, lYang et all 
d2012h and IMoster et al.1 |2012j) have pointed out the importance 
of using the infall redshif{j in HOD and ABM models. In our ap- 



1 The infall redshift is denned as the last time the satellite galaxy was the 
central object within its FOF group, see Eq.[T]below. 
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proach, the stellar mass does not depend directly on the infall red- 
shift. However, for a given satellite subhalo mass, the host halo 
m ass correlates well wi th the infall redshift. Moreover, according 
to lNeistein et al.l J201 lal) the combination of both subhalo and halo 
masses reproduces the clustering of SAM galaxies better than using 
the subhalo mass and infall redshift. 

This work includes the following parts. In section [2] we de- 
scribe our methodology and explain how we model the SMF, CF, 
and WL. In section[3]we show the results for the mass relation be- 
tween haloes and galaxies, and explore the contribution from dif- 
ferent constraints. Lastly, we summarize our main findings and dis- 
cuss them in section [4] Throughout the paper we write 'Log' to 
designate 'Logio'. 



2 HASH MODELS 

In this section we provide details regarding the method and models 
used in this paper. Our approach assigns a stellar mass to a set of 
subhaloes from a large iV-body simulation. The set of subhaloes, 
their mass and location are described in the next three subsections. 
In the rest of the subsections we explain how the SMF, CF, and WL 
are derived, and how we define specific models within our formal- 
ism. For a few example models that fit both the SMF and CF, the 
reader is referred to paper-/. 



2.1 The cosmological simulation 

We use info rmation on haloes an d subhaloes from the Millennium 
simulation dSpringel et aT1l2005l) . This simulation follows 2160 
dark-matter particles (each with mass of 8.6 x 10 s /i _1 M ) within 
a box of length 500 /i _1 Mpc. It is based on a cosmological model 
with parameters (O m ,fi a, (Ts, h) = (0.25, 0.75, 0.9, 0.73), and 
includes 63 output snapsho ts spaced by w 250 Myr. The FOF al- 
gorithm jPavis et"aL l fl985T) was used to identify h aloes, and was 
then used as an input for the SUBFIND algorithm dSpringel et al.l 
1200 lh to identify subhaloes. The merger trees used here are those 
based on subhaloes, as described in lSpringel et al.l d2005h . 



2.2 Haloes and subhaloes 

The subhalo mass, Mh, is defined here to be the mass of all the par- 
ticles inside a subhalo, as identified by SUBFIND. Throughout this 
paper we use the infall mass of subhaloes, Afi n f a n, which should 
correlate better with the stellar mass of galaxies, 



if central within its FOF group 



M infa ll = 



(1) 



yM hiP (zinia,n) otherwise 

Here Zi n f a ii is the lowest redshift at which the main progenitoi0 of 
the subhalo Mh was the most massive within its FOF group, and 
Mh, P is the main progenitor mass at this redshift. 

In addition to Minfaii we will use the halo mass M200, which 
is defined as the mass within the radius where the halo has an over- 
density of 200 times the critical density of the simulation. The halo 
mass for satellite galaxies within a group is defined to be the halo 
mass of the central object of that group. In general, M200 includes 

2 Main-progenitor histories are derived by following back in time the most 
massive progenitor in each merger event. 



mass not only from the central subhalo within a group, but also 
from all (or some) of its satellite subhaloes. 

At each redshift we construct a catalog of subhaloes, including 
all subhaloes that are identified at this redshift, and additional pop- 
ulation of unresolved subhaloes (see below). The catalog includes 
the following information: Mi n f a n and M200 (for fixing the stellar 
mass); the position (for computing the CF, i)2,5t ; and the projected 
dark-matter density profile (for computing the WL signal, §2.61 . 
Unless stated otherwise, we assume that all observations are fixed 
at z — 0, and compute the subhalo catalog at z = 0. 



2.3 Satellite subhaloes 

We divide the population of subhaloes into three types, 

• central subhaloes: most massive subhaloes within their 
FOF group. 

• satellite subhaloes: all subhaloes except central subhaloes. 

• unresolved subhaloes: subhaloes that were last identified at 
higher redshift, and are added according to the dynamical friction 
formula (Eq.|2]below). All the unresolved subhaloes are also satel- 
lite subhaloes. 

Within the iV-body simulation, satellite subhaloes lose their 
mass while falling into a bigger subhalo, and thus might fall below 
the resolution limit used by the SUBFIND algorithm (20 particles). 
However, the galaxies that should reside inside these satellite sub- 
haloes might live longer, as they are more dense, and thus less vul- 
nerable to stripping. While constructing the merger trees, subhaloes 
that fall below the resolution limit are considered to be merged into 
the central object. Consequently, this effect can modify the abun- 
dance of subhaloes even at two orders of magnitude above the min- 
imum subhalo mass resolved by the simulation (e.g. lNeistein et all 
l2011ah . 

In order to take this effect into account, we add to our sub- 
halo catalog at z — some subhaloes that have merged with larger 
subhaloes at a higher redshift. This is done by modeling the time it 
takes a galaxy to fall into the central galaxy by dynamical friction. 
At the last time the satellite subhalo is resolved (to) we compute its 
distance from the central subhalo (d sa t), and estimate the dynami- 
cal friction time using the Chandrasekhar formula, 



tdf 



<*df 



1.17K4 2 at 



GM h , 2 m(l + M M /A4, 2 ) 



(2) 



Here Mh,i is the mass of the central subhalo, V v is its virial veloc- 
ity, and Mfj.,2 is the mass of the satellite subhalo. Once a satellite 
subhalo falls together with its central subhalo into a larger group, 
we update tdf for both objects according to the new central sub- 
halo. All subhaloes for which tdf + to is larger than the cosmic 
time at z = are added to the catalog of subhaloes at z = 0. When 
using the total mass of t he satellite galaxy instead of Mh.i, Qdf 
was found to equal ~ 2 1 Colpi et al J[l999l: iBovl an-Kolchin et al .1 
120081 : 1 Jiang etalj|2008l ; iMo et alj|20icb . Since we only use the 
dark-matter subhalo mass here, and we would like to be careful 
when adopting limiting assumptions within our models, we allow 
Qdf to vary between 0.1 and 10. 

In order to model properly the location of unresolved sub- 
haloes we use two options. First, we use the location of the most 
bound particle of the last identified subhalo (as was done by e.g. 
ICroton et alj2006h . This particle is followed over time even though 
its host subhalo does not exist anymore. The particle location at 
z = is taken directly from the iV-body simulation. As a second 
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option we derive an analytical model for the infall of a subhalo in- 
side a central potential, 

d=d S at(l~T P ) 1/9 . (3) 

Here r = (t(z = 0) — io)/idf is the fraction of time spent out 
of all the estimated dynamical friction time until z — 0, d is the 
distance we adopt at z — from the central subhalo, and p, q are 
constants. The three dimensional location of the subhalo is using 
the assumption of radial infall. For a full derivation of this model, 
the reader is referred to paper-/. When searching for the models 
that fit observations we always try both options: the location of the 
most bound particle, and the analytical model above. 

2.4 Modeling the stellar mass function 

In order to model the stellar mass function of galaxies we first com- 
pute the mass function of subhaloes as found in the Millennium 
simulation and divided into central and satellite components: 

dA c 



=(M 5nfeU ) = 



T/dlogM infall 
1 



(4) 



d 2 N 3 



#*(MnfaU, M 200 ) = ^"Tj " . jTT • (5) 

V dlog M in& ii dlog M 200 

Here V is the volume of the simulation box, and N c , N s are the 
numbers of central and satellite subhaloes respectively. For central 
subhaloes we need to compute the mass function as a function of 
A/i n f a ii only, because there is no obvious physical reason to assume 
a dependence on both Minfaii and M200 ( on average, the value of 
A/200 for central subhaloes is smaller by only ~ 0.08 dex than 
A/mfaii, with an RMS scatter of ~ 0.06 dex). The mass function 
for satellite subhaloes is saved as a function of both M^tali and 
M200, allowing us to model the stellar mass as a function of these 
two different variables. 

Next we would like to use the mass function of subhaloes to 
predict the stellar mass function of galaxies. When comparing to 
the observed SMF, we need to compute the number of galaxies that 
corresponds to a given range in m 4 . Since the stellar mass is a func- 
tion of the subhalo and halo masses, = m* (Minfaii; A/200), 
the number of galaxies within a given range in m* is the same 
as the number of subhaloes within the corresponding area in the 
(A/infaii, A/200) plane. A simple way to parametrize regions within 
the (Minfaii, A/200) plane is by using the boundaries of each region. 
We define to be the region boundaries for central subhaloes, and 
U" the boundaries for satellite subhaloes. To conclude, The num- 
ber of galaxies within a given range of m* is simply the integral of 



N = V 



r u Ui 

L 



■ dlog A/infall + 



(6) 



X 



VI (j>s dlog Minfaii dlog A/200 ■ 



2.5 Modeling the auto-correlation function 

The usual way to compute the CF is to first populate a specific list 
of subhaloes from the simulation with galaxies, and only then to 
compute their CF. Since the CF is based on counting the number 
of pairs of subhaloes/galaxies, each such computation is very de- 
manding in terms of computer resources. In our HASH approach 
we first compute the underlying number of pairs of subhaloes at all 
individual masses. Only at a later step we integrate the pair numbers 



in order to compute the CF of galaxies. This method is extremely 
efficient when exploring many models, allowing us to compute the 
CF for ~ 10 7 different models. 

We use the entire set of central subhaloes from the Millennium 
simulation and count the number of pairs into tp C c- 



[fall, r ) — 



1 d 3 A^ c 
V 2 dlog A/j^^jj dlog M 2 niall dlog 



•(7) 



Here A/4 lfall , M 2 n{all are the infall mass of the first and second sub- 
haloes in the pair, and r is the distance between these subhaloes 
within the x-y plane, taking into account the periodic boundary 
conditions (we use the x-y plane in order to compute the projected 
auto-correlation function, as is described below). In practice, we di- 
vide the range in LogA/i n f a ii and Log r into bins of 0.1 dex (except 
for A/i n f a ii < 10 12 Me for which we use bins of 0.02 dex), and 
save tfj cc as a multi-dimensional histogram. 

In a similar way we count the number of pairs for central- 
satellite and satellite-satellite subhaloes, 



^(-Mifaii, Af£fau, A/ 2 2 o, r) . 

VMMnfall,A/infall > M 2( )0, M 200 ,r) 



(8) 
(9) 



Note that for satellite subhaloes the number of pairs is saved as a 
function of both Minfaii and A/200. This is done in order to properly 
model the dependence of stellar mass on A/200 • 

The observed CFs used in this work are based on the auto- 
correlation function within bins of stellar mass (hereafter 'do- 
mains'). We therefore need to compute the number of pairs N p for 
all galaxies within a range of stellar masses. Similarly to Eq.|6]we 
need to integrate the different tp components over the correspond- 
ing region within the Minfaii — A/200 plane: 



diV P (r) _ y2 
dlogr 



V 



ipccd logA/ infa n + 



Ipss d 2 log Minfai] d 2 log A/200 + 



Tpcs d log A/i n f a ll dlog A/ 2 , 



(10) 



(11) 



The projected two-point auto-correlation function, w p (r), is 
then defined as the deviation in the number of pairs from the aver- 
age value per volume: 



w p (r) = 



L 2 N p (r) 
N 2 A(r) 



(12) 



Here A(r) is the 2-dimensional area covered by the radial bin r, L 
is the size of the simulation box in /i _1 Mpc, and is the number 
of subhaloes (see Eq.[6]l. The factor N 2 A(r)/L 2 corresponds to the 
number of pair galaxies with a separation r, in case of a uniform 
random distribution of objects inside a 2-dimensional box of size L 
(note that we use the convention of double counting each pair when 
computing N p , so A^ 2 corresponds to twice the number of unique 
pairs). 



2.6 Modeling weak lensing 

In order to compute the WL signal of our HASH models, we first 
compute the projected two-dimensional density profile for each 
subhalo within the Millennium simulation. This is done by pro- 
jecting the entire database of particles into the X-Y plane, so the 
projected density profile of each subhalo includes integration on 
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Figure 1. The lensing profiles for different subhalo masses. Symbols show 
the average lensing profiles for central and satellite subhaloes within the 
Millennium simulation with Mj n f a n and M200 as indicated (profiles of 
central subhaloes are those where the given mass is only iWinfall)- F° r com- 
parison we show in so lid and dashed lines the same quantities for central 
subhaloes, taken from lHavashi & White! J2008L note that for central sub- 
haloes A/200 is smaller than Afj n f a n by ~0.08 dex). For satellite subhaloes 
we use here ayf = 3, and location of unresolved subhaloes following their 
most bound particle. 

the full simulation box along the Z-axijf]- These profiles are then 
averaged out into bins of subhalo and halo masses, resulting in: 

E (M infaU , r) , E s (M infall , M wo ,r) . (13) 

Here r is the two-dimensional distance to the centre of the subhalo, 
and the subscripts s/c mark satellite and central subhaloes respec- 
tively. We save E SjC as a multi dimensional histogram, with bins of 
0.1 dex for each variable (except for Mi n f a n < 10 12 , where the bin 
size is 0.02 dex), as was done in the previous subsections. 

The observed weak lensing signal is a transformation of the 
projected density profile: 

A£ c (M infall ,r) = (14) 

2 r 

— / riE c (M infaU ,ri)dri-E c (M tafall ,r), 
r Jo 

where the first term on the right hand side is the average profile 
within r. Note that AE C is a linear function of E c , and is zero for 
S e that are constant in the variable r (the mass sheet degeneracy, 
iGorenstein et alj[l988l) . Consequently, integrating the density pro- 
file of each subhalo along the Z-axis of the entire simulation box 
does not affect the value of AE C . For satellite subhaloes, the deriva- 
tion of AE S uses the same identity as above, but with E s . We refer 
to AE as the 'lensing signal' or just 'WL signal' throughout the 
paper. 

As was done above (Eq. |6), in order to compute the model 
prediction for a range of m* we integrate AE C and AE S over the 
values of Mi n f a ii and M200 that correspond to each domain. The 
integration uses the domain boundaries U£ and Uf exactly like in 
Eq.[6] with the only difference in using AE instead of <j>. 

In AppendixlAlwe list the specific values of AE SjC for a large 
sample of profiles, allowing users to easily apply this data to their 

3 We use a grid in the X-Y plane with a cell size of 10 /i -1 Kpc. 
© 2012 RAS, MNRAS 000. mfE] 



own use. In Fig. [TJ we show a few examples of profiles for both 
satellite and central subhaloes. It can be seen that the profiles of 
satellite subhaloes are more flat than for centrals, but they coincide 
with the profile of the central object at r > 1 Mpc. For central 
subhaloes, a factor of ~ 100 in the subhalo mass results in a factor 
of ~ 10 enhancement in the lensing signal. Also, more massive 
subhaloes have more extended lensing profiles. 

We also show in Fig. [Tj that the profiles taken from 
lHavashi & White! J2008h a gree well with ours. These authors have 
computed the same quantity as we do here based on the same sim- 
ulation, but with a different technique, and only for central sub- 
haloes . Note that for central su bhaloes we use the Mjnfati mass, 
while lHavashi & White! J2008h use the M200 mass (on average 
M200 is smaller than Mtnfoii in ~ 0.08 dex for the same central 
objects). 

When constructing E s directly from the simulation we only 
use the location of unresolved subhaloes following their most 
bound particle. However, we have many models in which the lo- 
cation of these subhaloes is set by an analytic recipe (Eq. [3j. In 
order to simplify the process of estimating E s for such models, we 
do not compute the density profiles directly from the simulation 
(using the new location) but rather compute the density profile an- 
alytically, using the profile of the central object within the same 
group. This is done by computing the overlapping area between 
each radial shell of the satellite subhalo and the shells of the central 
object. Each shell of the new profile is therefore a combination of 
various shells from the central object, each with a different weight 
factor, according to the size of the overlapping area. As the central 
profile used here is saved in radial bins, this method is accurate only 
for spherical density profiles, so it might introduce some errors in 
our case. 

In order to test the method above, we compute analytically the 
profiles of unresolved subhaloes, with the target location being the 
same as is obtained by their most bound particle. We then check if 
the profiles extracted from the simulation do agree with our analytic 
derivation. In Fig.|2]we show the results of this test. It is evident that 
this correction is not entirely accurate, giving rise to an offset of up 
to a factor of 2 at the smallest scale. However, since the effect of 
unresolved subhaloes is rather weak on our results (see below), the 
low accuracy shown here should not modify our results strongly. 

The inconsistency of our approximation above means that 
satellite subhaloes are located in regions of high density, also in 
comparison to the average location within the host halo, at the same 
radial distant from the centre. Consequently, when using a fully an- 
alytical model, one needs to assume the locations of satellite sub- 
haloes, as well as their excess density. In our approach, we neglect 
the excess density only for unresolved subhaloes (which are a small 
fraction of all the satellite subhaloes), and only in cases where the 
location is not set by their most bound particle. 

2.7 The definition of a model 

In this study we use five stellar mass bins (termed 'domains') start- 
ing at logM0 = 9.27 and ending at 11.77, each with size of 0.5 
dex. These domains were chosen because they match those of the 
observed CF used here. To simplify our formalism, we will use the 
same domains to interpret the observations of WL and SMF. 

As discussed above, in order to make predictions based on our 
models, we need to make assumptions on how each domain in stel- 
lar mass is related to a region within the (Mi n f a n , M200 ) plane. Our 
assumptions on these regions and their corresponding boundaries 
Ui are as follows. 
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Figure 2. The average lensing profile for unresolved subhaloes of mass 
log iW; n f a ii within the range [11, 11.5], and log M200 between 14 and 15. 
Blue solid line shows the actual profile measured from the Millennium sim- 
ulation, based on the location of the most-bound particle. The red clotted- 
dashed line shows the average profile of our estimated correction, in which 
each profile of a satellite subhalo is being computed from the central object 
profile, using the distance of the satellite. Black dashed line shows the av- 
erage profile of the central object within the same FOF group for all these 
subhaloes. 



Uf is only limiting the domain boundaries. The combination of 
different domains, and the functional dependencies within each do- 
main are free from any parametrization. To emphasize this point, 
any function of the type: 



= / [a(d) log Minfeu + 6(d) log M 20 o] 



(17) 



is possible using our parametrization, with a, b depending on the 
specific domain d (i.e. a & b can depend slowly on Mj n f a n as well), 
and f(x) is completely free. This is different from ABM studies 
that assume a fix functional shape for / and zero b. 

Using domains in stellar mass, and corresponding regions in 
the Minfaii — M200 plane saves a large amount of computation time, 
since we perform our analysis separately on each domain, reducing 
the number of free parameters by a factor of roughly five. During 
the post-processing phase the code accepts only models for which 
all the five domains have boundaries that coincide with their neigh- 
bour domains (for example, the upper boundary of the first domain 
is the same as the lower boundary of the second domain). 

Even though our parametrization introduces a large degree of 
freedom, our current approach does not include the effect of ran- 
dom scatter when assigning to subhaloes. Many works have 
claimed that a log-normal scatter with a ~ 0.15 — 0.2 is consistent 
with the observations. In our case, we choose to deconvolve the 
scatter from the observational results, a method that is equiv alent 
to adding a scatter to the model (see e.g. lBehroozi et alj2010h . We 
will use this method below, for checking the effect of scatter on the 
SMF and WL fits. 



For central subhaloes (where we assume that m* depends only 
on Mjnfa.11), each domain is defined as a range in Minfaii • This 
means that a range in stellar mass for central galaxies is equiva- 
lent to a range in Minfaii for central subhaloes. In other words, this 
parametrization simply assumes that the shape of the m* — Mi n f a n 
curve is determined by a sample of a few points, allowing its shape 
to freely depend on the data. In order to define a model we therefore 
need to fix one free constant for each boundary: 



Mm 



(15) 



Since we do not adopt prior limitations on M^ faI1 { this description 
does not restrict the models in any sense. 

For satellite subhaloes, we assume that each domain in stellar 
mass corresponds to a region of subhaloes that are located between 
two boundaries within the Minfaii — M200 plane: 



log U t s = log M^feu = + Si log M 2 i 



(16) 



where M^ lfall s and Si are free constants, M200 is a variable and U£ 
corresponds to the value of Minfaii (in analogy to the equation of a 
straight line, y = a + bx). This means that the boundaries of each 
domain are linear lines within the log Minfaii — log M200 plane. 
The values for Mj S nfall j are restricted to follow our basic sampling 
bins. For Si we write Si = tan 9i and sample 8i in steps of 6.75 
degrees, over the range [-90, 45]. Note that ABM models assume 
that Si are all zero, and Mnf a ii,i = •^infaii,i- 

The parametrization we adopt for satellite subhaloes is us- 
ing some prior assumptions that might restrict all the possible 
functional dependencies of m* (Minfaii, Afeoo). A more general 
model would allow each boundary to be a free contour within the 
Minf a ii — M200 plane. Such models, however, are hard to handle 
as the number of possibilities for each domain would be huge (in 
our simplistic case, and using the sampling density of our param- 
eter space, there are already ~ 10 10 possible boundaries per each 
domain). Nonetheless, our assumption on the functional shape of 



3 RESULTS 

3.1 Search strategy and fitting criteria 

In this section we conduct a detailed search within all the possi- 
ble HASH models, to see which model can fit a given constraint. 
Our search goes over all the possible boundaries Uf and U°, taking 
into account all the options for M^aii u -^infali i: Si as defined in 
Eqs.[T5]and[T6] In addition, we have tested values of ayf between 
0.1 and 10, and a few options for the location of unresolved sub- 
haloes (using either the location of the most bound particle of each 
subhalo, or various values of p and q from Eq. [3] as summarized 
in table Al of paper-/). In total we have 21 parameters for each 
model: 3 parameters that define each boundary times 6 different 
boundaries, and 3 general parameters that are assumed to be fixed 
over all domains: ayf, p and q. Overall, our search algorithm scans 
~ 10 11 models for each domain. 

When matching the SMF, we require an accuracy that is bet- 
ter than 20 per cent in reproducing the number of galaxies within 
each domain (Eq.[6]l. Once the model reproduces the total number 
of galaxies within each domain (of size 0.5 dex in m*) it is then al- 
ways possible to distribute the galaxies properly within this domain 
to fully reproduce the shape of the SMF. This is doable because 
our pa rametrization is not bein g limited within each domain. In ad- 
dition, iLeauthaud et alj d201ll ) have shown that the SMF errors at 
different masses are highly correlated, supporting the concept of 
fitting the SMF with only a few data points. 

For the CF and WL we estimate the quality of the models us- 
ing a simple RMS test: 



RMS, 



1 n 

l -y^ogwp,* -logwp.i) 2 , 

\ i=i 



(18) 



where w p .i is the observed CF value at the i-th point, w p ,i is the 
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Figure 3. The projected lensing signal derived for the models that fit WL only. Blue solid lines represent the maximum and minimum of all the models that 
fit the observed WL values to a level of [0.35 0.2 0.1 0.1 0.1] dex RMS for each domain respectively. The observational reference from M06 binned into our 
stellar mass domains are shown as symbols, with error bars that correspond to 1 standard deviation. 



model prediction, and the index i goes over the points measured for 
one domain. A similar definition is used for estimating the quality 
of the WL signal. 

We note that when using log values of CF and WL, the obser- 
vational errors are roughly constant within each domain, making 
our RMS test similar to the popular \ test. The benefit in using 
an RMS test is that our results do not depend on the specific errors 
given by the observational studies used here, allowing for an easy 
adaption for different data sets, and simple comparison of errors in 
WL versus CF. In addition, the RMS criterion is not sensitive to the 
number of data points for each constraint, unlike the \ 2 test that 
includes all data points together. Lastly, systematic uncertainties in 
our methodology (e.g. the specific cosmological parameters used 
by the simulation; the level of accuracy inherent in using two pa- 
rameters to fix m*) are not negligible, and it is not easy to estimate 
the full error bar for a given data point. 

It is well known that vario us errors within t he ob served data 
sets are correlated. For example. iLeauthaud et aljfeoill) have stud- 
ied these correlations in detail for the observations of SMF, CF, 
and WL. They show that CF and WL values at different neighbour- 
ing scales are correlated, as well as SMF errors for neighbouring 
masses. In addition, correlations might exist between different ob- 
servations. For example, the WL signal an d CF will both b e mod- 
ified in a similar way by the value of f2 m dYoo et ai]l2006t ). Since 
we do not have the full covariance matrix that describes correlations 
within all the data points, we choose to keep our fitting procedure 
simple and transparent, and use the RMS criterion above. 



To summarize, we list below all the criteria that are used to 
select models: 

• The number density of galaxies within each domain agrees 
with the observed value at a level of 20 per cent. 

• The RMS estimate for the CF is computed by taking the ob- 
served points at 0.03 < r < 30 Mpc h^ 1 , spaced by 0.2 dex in 
log r. We use two different criteria, of RMS=0. 1, 0.2 dex. 

• The RMS estimate for the WL is computed by taking all the 
observed points at 0.03 < r < 2 Mpc h^ 1 . Here we use two 
different criteria, an RMS level of [0.35 0.2 0.1 0.1 0.1] dex per 
each domain, or an RMS value of 0. 1 dex for all domains. 

• The number of galaxies more massive than the most massive 
domain (i.e. galaxies more massive than 10 11 ' 77 Mq) can deviate 
by no more than 20 per cent, with respect to the nominal value 
plus/minus an additional poisson error. This larger error is adopted 
because there are usually only a few tens of galaxies within this 
range. We note that even though galaxies at this mass range are not 
part of our five domains, we in practice fit the SMF at this mass 
range, similarly to what is done in ABM. 

• Models for which the fraction m*/Mi n f a n is bigger than the 
universal fraction of 0.17 are rejected. 

3.2 Models that fit weak lensing (WL) 

We first show how the WL can constrain the mass relation between 
haloes and galaxies when it is used as the only constraint for the 
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Figure 4. The mass relation between haloes and galaxies using constraints 
from WL only. Thick lines represent the maximum and minimum of all the 
models that fit the observed WL values to a level of [0.35 0.2 0.1 0.1 0.1] 
dex RMS for each domain. Results from M06 are plotted in symbols, and 
are shown only for central galaxies. Squares and circles refer to the mean 
value of early and late type galaxies respectively, with error bars that reflect 
95 per cent confidence level. The fraction of late-type galaxies out of the 
full sample is 0.74, 0.60, 0.46, 0.32, 0.20, 0.11, 0.05 (ordered in increasing 
m*). For reference we plot in dotted line the A4i n f a ii - m* relation using 
both satellite and central galaxies as one population , with no dependenc e on 
M20O: an d constrained by matching the SMF from|Li & White 1 2 0091) . 

models. The observation s used for the search are the data from 
Mandelbaum et al hereafter M06), re-binned to match the 

domains in stellar mass used for measuring the CF (these domains 
typically include two mass domains from M06jf|. The observed 
data from M06 after re-binning can be seen in Fig. [3] Since the 
two lowest mass domains are rather noisy, we allow their fit to de- 
viate by up to 0.35 and 0.2 dex RMS respectively. For the other 
domains, we demand an RMS fit of 0. 1 dex, which roughly agrees 
with an error of one standard deviation (the plotted error bars). 

Technically, in addition to the WL constraint, we have re- 
quired a factor o f 50 accuracy in matching the SMF from 
iLi & White! d2009T) . Unless we apply such a limit, the number of 
models gets too large for us to analyze. Nonetheless, this constraint 
is rather weak, and should not affect the model results too strongly. 

The results of our search in terms of the MR are plotted in 
Fig. [4] where we plot the maximum and minimum Mi n f a n for a 
given m*, using all the accepted models. Note that for central sub- 
haloes each HASH model predicts a unique MR, with no scatter. 
For satellite subhaloes, each HASH model might have a range of 
A^infaii per a given m*, depending on the mass of the host halo, 
M2oo- We therefore compute the median value of Minfaii at each 
domain boundary Ui (m*=9.27, 9.77, 11.77), and then plot 
in Fig. [4] the maximum and minimum for all these median val- 
ues. The set of models that are able to match the first domain 
(9.27 < log(m*/M0) < 9.77) is too large for us to analyze, 
and is therefore not shown here. 

Regarding central galaxies, our HASH models predict a sig- 



4 When re-binning the data from M06 we take into account the number of 
lenses within each original bin. We have also tested our search algorithm 
using the original bins from M06, finding no significant difference. 



nificant level of freedom, of more than one order of magnitude un- 
certainty in stellar mass, for a given subhalo mass. Consequently, 
the WL signal is not useful as a single constraint. We will show in 
the next section that this is also true when combining constraints 
from WL and SMF. For satellite galaxies, our search shows even 
larger uncertainties, probably because the profiles of satellite sub- 
haloes have lower values than for centrals, and less dynamical 
range (see Fig. [T). Another hint for this is obtained from compar- 
ing the range of models above, against the range of models obtained 
where the locations of unresolved subhaloes are fixed at the posi- 
tion of their most bound particles. This reduction in the parameter 
space shows a small effect on the mass relation, hinting that unre- 
solved subhaloes are not a key ingredient in fixing the weak lensing 
signal. 

Our results from Fig. [4] are consistent with the findings of 
M06 (shown for central galaxies). The fraction of massive late type 
galaxies (plotted in circles) is lower than ~ 10 per cent, and has a 
negligible effect on the lensing signal. As a result, the large error 
bars for these galaxies do not point to a discrepancy in compari- 
son to our results. It is interesting that our range of models is much 
larger than that of M06. It seems that the methodology used here 
introduces more freedom in the models. Another option is that our 
method of directly scanning the parameter space has some benefits 
over the Monte Carlo Markov chain approach used by M06. How- 
ever, our fitting criteria are somewhat different from those adopted 
by M06, making it hard to compare the two approaches directly. 
Note that the cosmological model assumed here is similar to the 
one assumed in 

In order to demonstrate the range of possible models that are 
able to match the WL, we show in Fig.|5]two different models that 
fit the WL as discussed above. As can be seen from this plot, dif- 
ferent models might have very different contributions from central 
and satellite galaxies. In the 2nd domain, model A is dominated by 
centrals at small scales, while model B is dominated by satellites. 
At the most massive domain, the difference between the models is 
maximal, and each model is governed by only one population of 
galaxies (central/satellite). 

3.3 Models that fit the SMF and WL 

The observed WL signal provides a relatively poor constraint on 
the MR. It is therefore interesting to see if adding the SMF would 
improve the result s. Since we are going to use stellar masses from 
lLi & White! ( [20091 hereafter Li09) in the next section, we first apply 
their SMF constraint together with the WL of M06. Interestingly, 
even though the range of models that fit the WL is relatively large, 
we do not find any model that matches the SMF constraint as well. 

The difference in stellar masses between Li09 and 
iKauffmann eTaU d2003l . used by M06 for deriving the WL 
signal) was examined in the Appendix of Li09. It is shown there 
that the two mass estimates deviate by roughly 0.1 dex, where 
Li09 tend to assign lower masses for the same objects. This might 
contribute significantly to the discrepancy between these two data 
sets found here. In order to bypass this issue we construct a model 
search in which we fit the WL from M06 toge ther with the SMF 
derived directly from [Kauffmann et alj d2003l . we use the SMF 
version presented in Li09). The results of this search are shown 
in Fig. [6] The set of accepted models looks very narrow, and the 

5 M06 have assumed the following set of parameters: (fi m , ag) = 
(0.3, 0.7, 0.9) in comparison to (0.25, 0.75, 0.9) here. 
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Figure 5. The projected lensing signal for two example models that fit WL (as summarized in Figs.[3]and[4). Thick solid lines are the total signals per each 
domain, dashed and dotted-dashed lines show the contribution to the profiles from satellites and central subhaloes respectively. Each model is represented by 
a different line color. For the most massive domain, each model includes only one population of subhaloes (i.e. only central or only satellite subhaloes). For 
clarity, we do not show solid lines for this domain. 



models have high fraction of satellite galaxies at the massive end 
(for a given value of m*, the value of A/i n f a n for satellite galaxies 
is much larger than for centrals). 

As will be discussed below, when using a mock WL signal, 
our results indicate to a much larger range of models, in compari- 
son to the models obtained when using WL from M06. This issue 
raises the concern that our models suffer from a restrictive or wrong 
assumptions. We have made the following tests in order to study the 
possible effect of our assumptions: 

• The mean redshift of lens galaxies from M06 goes up with 
stellar mass, reaching z = 0.19 at the massive end. On the other 
hand, our analysis assumes that the observations are all referring to 
z = 0. We have tested the importance of this effect by using lensing 
profiles of subhaloes taken at the proper redshift for each domain 
(while using the same SMF from low redshift). Results for the MR 
in this case are similar to Fig. [6] However, a full exploration of this 
issue will require to allow the boundaries of neighbouring domains 
to differ, because the MR might evolve with redshift. Taking this 
effect into account deserves a further analysis that cannot be done 
here. 

• We have tested the range of models when using the original 
bins from M06, based on either the stellar mass, or the luminosity 
(lum inosity based WL i s combined with the luminosity function 
from lBlanton et alj2003h . In both cases the range of models is sim- 
ilar to what we plot in Fig. [6] 

• We have tested the effect of adding a random scatter to the 



value of in-),. Assuming the scatter does not depend on Mi n f a ii, we 
have deconvolved the observed SMF and WL with a constant log 
normal scatter of 0.1, 0.15, and 0.2 dex. For each different scatter 
parameter, this results in a new SMF, and a new WL signal (we 
deconvolve the WL signal as a function of m* for each radial bin, 
taking into account the averaging of AE over the domain range 
in m*). These data sets can then be fitted by our models with no 
scatter. In order to check the effect of this on the consistency be- 
tween the SMF and WL, we have computed new domains in 
that will yield the same number density of galaxies as in the orig- 
inal domains, but using the new deconvolved SMF. For these new 
domains we also generate the WL signal, using the deconvolved 
signal above. Interestingly, even though both the SMF and WL 
change due to the addition of scatter, the new domains have very 
similar WL profiles to the original domains, with deviations that 
are smaller than our fitting criteria. Since our search procedure does 
not depend on the actual stellar mass value of each domain, it will 
produce the same solutions using the new domains above, as were 
found with the original domains. The only effect of scatter would 
thus be to shift the domains of m* to lower values, especially at 
the high mass end. This will have no effect on the range of models 
presented in Fig. [6] 

• The Millennium simulation used by our code assumes a 
higher v alue of as than the latest estimates (0.9 instead of 0.8, 
see e.g. iKomatsu eta7]|2009h . The effect of various cosmologi- 
cal models was examined closely by a recent set of papers using 
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Figure 6. The mass relation using constraints from both WL and SMF 
separated to central and satellite galaxies (upper and lower panels respec- 
tively). Fitting criteria for matching the WL signal are the same as in Fig. [4] 
The SMF is fitted to a level of 20 per cent for each domain. Lines corre- 
sp ond to models that are c onstrained by WL from M06, and SMF based 
on iKauffmann et alj <2003l) . The shaded regions correspond to models that 
fit the SMF from Li09, together with modified WL signal (the modified 
WL signal is chosen to agree with Li09, see section l3~4"t . The dotted line is 
plotted for reference, and is the same as in Fig. [4] 



the CLF approach l lvan den Bosch et a i] |20lllMoreetalJl2012l; 
Cacciato et all2012h . From these works (see especially iMore et al.l 
2012b it seems that the combination of SMF and WL is less sensi- 
tive to the cosmological model than the combination of SMF and 
CF. As we will show below, our model is able to match both the CF 
and SMF to a good accuracy. This might indicate that the cosmo- 
logical model assumed here is not the reason for the b ias we detect 
in the WL signal. However, since IMore et alj 12012b have used a 
different model than what we use here, this conclusion should be 
taken with a grain of salt. A clear test to this issue can only be 
done by using a different iV-body simulation with more accurate 
cosmological parameters. 

• We define a 'modified' WL signal, by adopting fiducial WL 
data points that are fully consistent with the observed CF & SMF 
constraint used below. Although this WL signal is arbitrary, it pro- 
vides a way to estimate how the WL constraint will affect our re- 



Figure 7. The mass relation for central galaxies, using various sets of ob- 
servational constraints. All models fit the SMF to a level of 20 per cent. 
The labels CF(0.1) and CF(0.2) refer to models that fit the CF to a level of 
0.1 and 0.2 dex RMS respectively. The label mod-WL corresponds to mod- 
els that match the modified WL signal to a level of [0.35 0.2 0.1 0.1 0.1] 
dex RMS, while mod-WL (tight) designate models that fit the modified WL 
signal to better than 0.1 dex RMS for all domains. All shaded regions cor- 
respond to the maximum and minimum M; n f a n values for each given m 4 . 
The solid line is given for reference and is the same as in Fig. [4] 



suits, in case it is fully consistent with the other observables. We 
explain how we derive the modified WL signal in the next section. 
In Fig.|6]we show the results of our model search using SMF based 
on Li09, together with the modified WL signal. For consistency, 
we use the same error estimate for the WL, as was done above 
when using the data from M06. This means that we demand an 
RMS accuracy of [0.35 0.2 0.1 0.1 0.1] dex for each domain. Note 
that the range of accep ted models when using M06 together with 
IKauffmann et all d2003T) is much more narrow than the range when 
using the modified WL signal together with Li09. 

We conclude that the discrepancy between the SMF and WL (which 
is also valid when using the CF below) might be due to the limi- 
tations of our study. However, it might also be that some observa- 
tional systematics contribute to this effect. 

From Fig.[6]it can be seen that using constraints from both the 
SMF and WL (the modified signal) is not so different from using 
only the WL (as seen in Fig. [4j. The range in Mi n f a n for a given 
stellar mass is still large for central galaxies, and even larger for 
satellite galaxies. We conclude that using both the SMF and WL as 
constraints on our models does not restrict the MR significantly. 



3.4 Models that fit the SMF, CF and WL 

In paper-/ we have made a search for all the HASH models that 
fit both the SMF and CF of galaxies. The range of models and the 
fitting criteria tested there, are exactly the same as being used here. 
The observed CFs are computed using the same stellar ma sses as in 
Li09, following the technique presented in iLi etal.|j2006f) . CFval- 
ues are given within five domains in log(m*/M0): [9.27, 9.77], 
[9.77 10.27], [10.27 10.77], [10.77 11.27], [11.27 11.77]. The CFs 
correspond to the projected two point auto-correlation function at 
scales ranging from 0.03 to 30 /i _1 Mpc (the two most massive do- 
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Figure 9. The projected leasing signal derived for models that fit the CF & SMF. Blue solid lines represent the maximum and minimum of all the models that 
fit the SMF to a level of 20 per cent, and the CF at the level of 0.1 dex RMS. Blue dashed-dotted lines show how the range of models changes, when the CF is 
constrained to fit the observations to a level of 0.2 dex RMS. The dots show the modified WL signal we adopt throughout the paper. These values are chosen 
to lie in the middle of the range of models that fit both the SMF and CF, at radial locations that are the same as in M06. The observational reference from M06 
are shown as symbols. 
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Figure 8. Same as Fig. [7] but for satellite galaxies. Here we first compute 
the median values of Af; n f a n for each m* value. The shaded regions de- 
scribe the maximum and minimum of the median values for all models. 

mains do not include data at small scales). Note that both the CF 
and WL signals go down to 0.03 h^ 1 Mpc. 

Computing the CF for each HASH model is much more time 



consuming then computing the WL signal, because we need to in- 
tegrate functions of higher dimensionality. As a result, we could 
not run a test in which only the CF is used as a constraint, like was 
done here with the WL signal. There are two searches being done, 
each using both SMF and CF as constraints, but with fitting the CF 
to an accuracy of 0.1 and 0.2 dex RMS. In Figs.[7]and[8]we show 
the MR for these two different criteria. In comparison to using the 
WL and SMF, we see that here the MR is more tight, especially for 
satellite galaxies. This however might be due to either the smaller 
errors used to fit the CF (note that the observed CF values are based 
on SDSS DR7, while M06 is based on DR4) or the larger range of 
scales probed by the CF. 

To summarize the results of fitting the CF and the SMF from 
paper-/, the range of models highly depends on the accuracy by 
which we fit the CF. Models that fit the CF to a level of 0.2 dex 
RMS show a large range of MR, especially for low mass cen- 
tral galaxies, reaching a factor of ~ 30 uncertainty in Mi nta n for 
log(m*/M ) = 9.27. This demonstrates that using the SMF as the 
only constraint does not yield tight MR within our formalism. De- 
manding higher accuracy in matching the CF (to a level of 0.1 dex 
RMS) improves significantly the MR for central galaxies, but has a 
relatively minor effect on the MR for satellite galaxies. Although an 
accuracy of 0. 1 dex RMS is larger than the quoted observed error 
bars (these are usually between 0.03 and 0.06 dex, depending on 
the domain and on the separation distance), we could not find mod- 
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Figure 10. The fraction of satellite galaxies out of all galaxies as a function 
of stellar mass. The maximum and minimum values for all the models that 
fit each set of constraints are shown in different line types as indicated. 
Labels are the same as in Fig. [7] 



els that match the CF to better than 0.08 dex. This probably means 
that systematics errors contribute significantly to the fit quality. For 
more details on the models that fit the CF and SMF the reader is 
referred to paper-/. 

In order to understand better the interplay between the CF and 
WL constraints we compute the predicted WL signal for all the 
models that fit the SMF and CF. The results of this test are shown in 
Fig. [9] For the three low mass domains, the WL signal predicted by 
our models looks similar to the observed data from M06. However, 
at the two highest mass domains, the models that fit the SMF & CF 
are actually inconsistent with the observations of WL from M06. 
Interestingly, the discrepancy still exist if we allow for a larger de- 
viation from the observed CF, at the level of 0.2 dex. All this agrees 
with the discrepancy between WL and the SMF discussed above. 
It should be noted that some part of the discrepancy seen in Fig. [9] 
especially at high stellar masses, is due to the different stellar mass 
estimates used in Li09 in comparison to M06. 

When combining constraints from SMF, CF and WL, the com- 
bination of the current data sets and our models is inconsistent. We 
have chosen to overcome this issue by modifying the WL data. In- 
stead of using the data from M06, we will now define a new mock 
WL data set, that agrees with the CF and SMF. At each domain we 
define the 'modified' WL signal by computing the average signal in 
log values, using the upper and lower limits of all the models that 
fit the SMF and the CF (we use a CF fit of 0. 1 dex for this purpose). 
This modified WL signal is computed at the same radial bins as the 
original M06 data, and is plotted in Fig. [9] 

In Figs.|7]and[8]we show the MR when using constraints from 
the SMF, CF, and WL. We use the SMF from Li09, CF as described 
above with an RMS fitting accuracy of 0.1 dex, and the modified 
WL signal with fitting accuracy of [0.35 0.2 0.1 0.1 0.1] dex RMS. 
This combination reflects the current accuracy of the various data 
sets. As can be seen from these figures, the range in MR is almost 
the same as the range when using only the SMF and CF. We con- 
clude that with the current uncertainties of the observed WL signal, 
it adds very little to our knowledge of the mass relation between 
galaxies and subhaloes. 

The error bars we assume for the WL signal are large, and are 



related to the fact that our reference data from M06 is based on DR4 
from SDSS. In order to estimate the potential of future WL mea- 
surements we now demand a tighter fit to the WL, at the level of 0. 1 
dex RMS for all domains. From Figs.[7]and[8]it seems that such an 
improved accuracy in the WL reference can significantly improve 
the accuracy in the MR for low mass galaxies. However, the MR for 
massive satellite galaxies remains uncertain even in this extreme 
case, with an order of magnitude freedom in the value of Mj n f a ii 
for lo fffm^/Mp)) = 11.77 (see the discussion in Ivan Uitert et al.l 
l201lL related to WL). 

An accuracy of 0.1 dex is 10 times higher than the value of 
0.35 dex we use with the current data. Assuming the errors in the 
WL signal are mostly statistical, this implies a factor of 100 in- 
crease in the number of galaxies in comparison to SDSS DR4 used 
here. Such a large number of galaxies is not expected to be ob- 
served in the near future, but might be doable with future surveys 
(e.g. EUCLID). 

Let us compare the power of the CF versus the WL in con- 
straining the models. From the solid lines in Fig. [9] (the most mas- 
sive domain) it is evident that for galaxies with m* > 5xlO n M Q , 
an error of 0. 1 dex in fitting the CF, corresponds to 0.05 dex range 
in the WL signal. This means that for these massive galaxies, us- 
ing a constraint of 0.1 dex for the CF (combined with the SMF) 
is equivalent to using a 0.05 dex error in the weak lensing signal. 
Consequently, in order for a new WL measurement to be power- 
ful, it needs to have much tighter error-bars than the existing CF 
measurements. 

In Tables [T] and [2] we summarize the MR values obtained for 
the different constraints used in this work. 

3.5 Satellite fractions 

In Fig.[l0]we show the fraction of satellite galaxies out of all galax- 
ies in the five different stellar mass domains used here. As was 
discussed above, models that are constrained by WL show a very 
large range of fractions, spanning values between 0.05 and 0.8. This 
is a much larger range than what has been obtained by previous 
studies (e.g. M06), again indicating that our formalism is proba- 
bly more general than other existing models. When using all con- 
straints (i.e. SMF, CF and the modified WL), the satellite fraction 
are more constrained, reaching a range of ~0.2 to ~0.4. This range 
is l arger than other estim ates, especially at the high mass end (see 
e.g. lCacciato et alj2012h . 



4 SUMMARY AND DISCUSSION 

In this work we have developed a formalism (termed 'HASH') to 
interpret various observations regarding the mass relation between 
haloes and galaxies (the 'MR'). Our approach is using a set of 
subhaloes from a large cosmological A-body simulation, assign- 
ing one galaxy to each subhalo. The stellar mass of galaxies is as- 
sumed to depend on both the subhalo and host halo masses (Afi n f a n 
and M200 respectively), allowing us to treat differently satellite 
and central galaxies. We devote a specific attention to satellite sub- 
haloes that are stripped within the A-body simulation, but might 
still correspond to observed galaxies. We allow some freedom in 
the number of these subhaloes via a model of dynamical friction, 
with a free constant scaling factor. The location of these unresolved 
subhaloes has a non-negligible degree of freedom as well. 

Our approach tries to adopt as much information from the sim- 
ulation as possible, while allowing freedom in quantities that are 
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Table 1. The possible range in the host subhalo mass M; n f a n, for central galaxies of a given m*, as derived from the various constraints used here (the same 
information is plotted in Figs.[4]|6]&|7). Column labels are the same as in Fig. [7] All columns other than the first one correspond to the subhalo mass, M; n f a n. 
All mass units are log Mq 





WL M 06 


SMF+WL mod 


SMF+CF0.2 


SMF+CFo.i 


SMF+CFo 


i+WL mod 


SMF+CF . 1+ WL mDd>tight 


9.27 






11.02- 12.50 


11.12- 11.96 


11.12 


- 11.78 


11.28- 11.34 


9.77 


10.98 - 12.50 


11.24- 12.26 


11.26- 12.60 


11.34- 11.98 


11.36 


- 11.96 


11.50- 11.58 


10.27 


11.54- 12.70 


11.62-12.70 


11.64- 12.70 


11.72- 12.20 


11.74 


- 12.10 


11.80- 11.94 


10.77 


11.98- 13.20 


12.40-13.40 


12.40- 12.80 


12.50- 12.70 


12.50 


- 12.70 


12.50 - 12.60 


11.27 


12.60 - 15.60 


13.70-15.60 


13.70- 15.60 


13.70- 13.80 


13.70 


- 13.80 


13.70- 13.80 


11.77 


14.40 - 15.60 


15.00-15.60 


15.00- 15.60 


15.00- 15.60 


15.00 


- 15.60 


15.00 - 15.60 



Table 2. Same as Table[JJ but for satellite galaxies. Note that here we show the median M; n f a n value for each m*, computed by going over all the possible 
values of M200- 



m* 


WL M 06 


SMF+WL mod 


SMF+CF0.2 


SMF+CFo. i 


SMF+CFo 


i+WL mod 


SMF+CFo . i +WL mod , t i ght 


9.27 






10.80- 11.70 


10.80- 11.47 


10.90 


- 11.47 


11.32- 11.45 


9.77 


10.84 - 12.24 


10.83- 14.30 


10.96- 11.98 


11.08- 11.75 


11.16 


- 11.75 


11.59- 11.70 


10.27 


10.98 - 12.79 


10.84- 12.42 


11.35 - 12.36 


11.42- 12.10 


11.56 


- 12.10 


11.92- 12.10 


10.77 


11.15 - 13.43 


11.06- 13.14 


12.03 - 12.98 


12.07- 12.81 


12.20 


- 12.81 


12.47- 12.81 


11.27 


11.66- 15.60 


11.66-15.60 


12.90- 14.10 


12.90 - 14.06 


12.90 


- 14.06 


12.90 - 14.06 


11.77 


11.66- 15.60 


13.40-15.60 


13.40- 15.60 


13.40- 15.60 


13.60 


- 15.60 


13.60- 15.60 



not accurately modeled. We thus hardly assume any prior on the 
shape of the MR. The only limitation built into our formalism as- 
sumes that the stellar mass of satellite subhaloes is a function of 
the linear combination of log Minfaii and log M200, and this only at 
six different values of stellar mass. Therefore, the MR for satellite 
and central galaxies vary significantly for different models within 
our formalism. In addition, we are able to scan systematically the 
full parameter space to a high resolution, consistent with our fit- 
ting criteria (~ 10 12 models), avoiding the uncertainties related to 
complicated search algorithms within the parameter space. 

The observational constr aints used here are: the stellar mass 
funct ion (SMF, we use both iKauffmann et alj|20ol lli & White] 
1 20091) ; the two-point auto correlation function (CF, we use a new 
version of lLi et all 120061 based on SDSS DR7); and a measure- 
ment of galaxy- galaxy weak gravitationa l lensing (WL, based on 
SDSS DR4 from lMandelbaum et all2006h . 

We claim that our formalism allows for more freedom in 
fitting the observed data sets than other methods. We explic- 
itly show that when matching the WL signal, our method pre- 
dicts a huge amount of freedom in the MR, reaching a factor 
of 100 uncertainty in Mi n f a ii for a stellar mass of ~ 1O 11 M0 
(see Fig. |4j and number-fraction of satellite galaxies that lie be- 
tween 0.05 and ~ 0.8 for all stellar masses (Fig. [10]!. These re- 
sult are very different than what was claimed in the past based on 
different approaches dGuzik & Selialdl2002t iHoekstra et alj|2004l : 



iMandelbaum et al.l2006l : lHoekstral2007l ; lvan Uitert et al.l201lh . al- 
though fitting criteria somewhat differ between different studies. 
Another example is that our models are hardly being constrained by 
the SMF, in contrast to halo occupation distribution (HOD) mod els 
for which the SMF plays a key role (e.g. iLeauthaud etlil]|20 1 ll) . 

Since all models used in this field are purely empirical, most 
of the assumptions being made are motivated by either more com- 
plex models, like hydrodynamical simulations, or by observational 
trends. Nonetheless, different assumptions might highly affect the 
generality of the model, and might induce significant limitations 
to the set of solutions found to match the data. Therefore, it might 
be that extending our formalism would further increase the set of 
models that can fit the data. On the other hand, other existing mod- 
els might suffer from a limiting set of assumptions as well. For 



example, HOD models assume that the number of satellite galaxies 
within a halo scales as a power law in the halo mass. This might 
not be related directly to the physical processes of galaxy forma- 
tion, and might not be valid in all the possible scenarios of galaxy 
formation. We argue that a much careful treatment of the assump- 
tions made by the models should be considered, in order to obtain 
the most general formalism possible. 

The HASH formalism used here is specifically different from 
most previous studies in the way the WL signal is being mod- 
eled (for more similar works see iTasitsiomi et alj |2004h iLi et alj 
120091) . Here we compute the density profile within the simulation 
around each subhalo. We then stack all profiles of the same sub- 
halo and halo masses (Mi n f a u, Af2oo)- This in principle provides 
a more accurate model than the usual, analytical HOD approach, 
at least for massive haloes that are well resolved within our sim- 
ulation. According to our analysis, the observed WL signal from 
IMandelbaum et alj 12006]) is not fully consistent with respect to the 
SMF based on the same stellar masses. We suspect that this is due to 
either the specific cosmological model assumed here, or to the non- 
uniform redshift of the data that cannot be modeled easily. Conse- 
quently, when combining constraints from SMF, CF, and WL, we 
use a 'modified', mock WL signal, that agrees with both the SMF 
and CF used here. 

In iNeistein et al. we have applied the constraints of 

both the SMF and CF, and showed that a large range of HASH 
models can match the data. The degeneracy in those models corre- 
sponds to ~ 0.8 dex freedom in the MR for low mass galaxies and 
for satellite galaxies of all masses. When we apply the additional 
WL constraint here, this range of models is hardly modified (see Ta- 
bles[JJand|2]for exact numbers). Based on the modified WL signal 
we claim that future WL observations with an uncertainty of ~ 0.1 
dex, would be powerful in constraining the mass relation of low 
mass galaxies. Future surveys like EUCLID are therefore crucial in 
reaching such high accuracies. For massive galaxies, it seems that 
constraints based on the CF are more powerful than those based on 
WL, for the same error estimates. 

There are various caveats within our approach. First, we use 
an 7Y-body simulation based on a non-accurate set of cosmological 
parameters and a limited volume (we use as = 0.9 in compari- 
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son to 0.8 predicted bv lKomatsu et alj|2009h . Second, our models 
do not allow for a random scatter in the MR (when computing the 
CF), which should naturally exist at some level. Third, assuming 
that the stellar mass of galaxies depends only on the mass of the 
host subhalo and halo is a significant simplification. From more 
complex models it is known that the large-scale environment might 
effect stellar ma s ses up to a few per cent (the 'assembly bias' effect 
iGao et al.|[2005l ; ICroton et al]|2007h. Also, the tim e a galaxy be- 
came a satellite might play a role T?ang et alj|2012h . Obviously, an 
accurate determination of the galaxy mass within a halo must take 
into account many properties of the merger history of that halo. 
Since all these effects are predicted to be at a level of a few per 
cent, our approach of fitting the data to a level of 0.1 dex might be 
secure enough. 

In terms of modeling the WL signal, our approach includes 
an additional set of assumptions. We assume here that the lensing 
signal is fixed only by the dark matter content. However, baryonic 
dynamics might change the dark-matter profiles of subhaloes, and 
their location within a group. In addition, t he mass of gas and star s 
can contribute directly to the lensing signal ^Leauthaud et al.l2012h . 
All these assumptions should be examined in the near future, as 
the observed data sets are getting more accurate, demanding more 
complex and general models. 
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APPENDIX A: DATABASE OF PROFILE SHAPES 

We provide the values of the WL profiles, AE C and AE S , averaged 
in bins of Mj n f a ii and M200 • The original tables have a bin size of 
0.1 dex in each variable (except for Mjnfaii < 10 12 for which the 
bin size is 0.02 dex), so we only write a sub-sample of bins here. 
For example, only half the values as a function of r are written. 
In addition, the number of subhaloes summed over all the lines 
within a table does not accumulate to the number of subhaloes in 
the simulation. 
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Table Al. The average WL profiles (AS C ) for central subhaloes, as a function of Af; n f a n. Units of A'/; n f a n are log Mg, profile units are /iMg pc 2 , and 
distances are given in log Mpc h . N is the number of objects used for deriving the average. 



Afinfall 


N 


log r = —1.5 


-1.3 


-1.1 


-0.9 


-0.7 


-0.5 


-0.3 


-0.1 


0.1 


0.3 


0.5 


0.7 


0.9 


11.0 
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Table A2. The average WL profiles (AS a ) for satellite subhaloes, as a function of Mj n f a n and M2oo- Here we use a model with a^f = 3, and location of 
unresolved subhaloes following their most bound particle. Mass units are log Mq, profile units are /iMq pc~ 2 , and distances are given in log Mpc h^ 1 . N 
is the number of objects used for deriving the average. 
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